arXiv:1505.02130vl [physics.chem-ph] 7 May 2015 


Dynamics and quantumness of excitation energy transfer through a complex quantum 

network 

M. Qin/ H. Z. Shen,i X. L. Zhao,i and X. X. Yi*^ 

^School of Physics and Optoelectronic Technology, 

Dalian University of Technology, Dalian 116024, China 
^Center for Quantum Sciences and School of Physics, 

Northeast Normal University, Changchun 130024, China 
(Dated: May 11, 2015) 

Understanding the mechanisms of efficient and robust energy transfer in organic systems provides 
us with new insights for the optimal design of artificial systems. In this paper, we explore the 
dynamics of excitation energy transfer (EET) through a complex quantum network by a toy model 
consisting of three sites coupled to environments. We study how the coherent evolution and the noise- 
induced decoherence work together to reach efficient EET and illustrate the role of the phase factor 
attached to the coupling constant in the EET. By comparing the differences between the Markovian 
and non-Markovian dynamics, we discuss the effect of environment and the spatial structure of 
system on the dynamics and the efficiency of EET. A intuitive picture is given to show how the 
exciton is transferred through the system. Employing the simple model, we show the robustness of 
EET efficiency under the influence of the environment and elucidate the important role of quantum 
coherence in EET. We go further to study the quantum feature of the EET dynamics by quantumness 
and show the importance of quantum coherence from a new respect. We calculate the energy current 
in the EET and its quantumness, results for different system parameters are presented and discussed. 

PACS numbers: 05.60.Gg, 03.65.Yz, 03.67.-a, 71.35.-y 


I. INTRODUCTION 


The processes of energy and information transfer in 
quantum networks play an important role in many areas 
of physics, chemistry and biology. They have been iden¬ 
tified as central to quantum mechanics since its early 
days [I|. Under realistic physical conditions, the cou¬ 
pling of quantum system to environment is usually un¬ 
avoidable, which leads to the deterioration of perfor¬ 
mance for fundamental mechanical processes in systems 
such as computing devices and biological organisms. Re¬ 
cent experiments, however, find evidence of long-lived 
quantum coherence in conditions that are often defined 
as hot and wet [3, Q in the Fenna-Matthews-Olson 
(FMO) protein complex of the green-sulfur bacterium 
Chlorobiumtepidum and in the reaction center of the 
purple bacterium Rhodobactersphaeraides 1M3. 

Inspired by these experimental observations, quan¬ 
tum coherence across multiple chromophoric sites has 
been suggested to play a significant role in achieving the 
remarkable efficient EET in EMO and other pigment- 
protein complexes [HHil- Since then, many works an¬ 
alyzing energy transport in these systems have been 
carried ou t 11414^ . such as by the quantum network 

model [3117 by the hierarchic equation [3 - l3| . by 

the renormalization group (3 . [3| , and by semiquantum 
theory [s^- Theoretical investigations found that dephas¬ 
ing does not always hinder the efficiency of EET com¬ 
pared with perfectly quantum coherent system, which is 
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contrary to the conventional intuition 31| . The interplay 
of the quantum coherent evolution and dephasing con¬ 
tributes to the highly effective EET in light harvesting 
complexes [3. f3|. Understanding the underlying mech¬ 
anism of such a process may assist us in designing novel 
nanofabricated structures for quantum transport and op¬ 
timized solar cells. 


Although there are many progresses made in this field, 
the interplay of the quantum coherent evolution and 
the relaxation process, however, has not been stud¬ 
ied in the Haken-Strobl model, in which the en¬ 
ergy relaxations are not included. In this pa¬ 
per, we introduce the relaxation process into the 
decoherence terms, and study the dynamics and 
the transfer efficiency in a complex network. We 
obtain several analytical solutions to illuminate how 
the coherent evolution and the noise-induced decoher¬ 
ence processes, especially the relaxation process, work 
together to reach efficient EET. 

Besides, the experimental evidence has shown that the 
space distribution of the pigments is very important for 
the EET 1^, Q . Thus it is essential to take phase fac¬ 
tors in the site-to-site couplings into account. We show 
how the spatial structure of the system, embodied in the 
site-to-site couplings [3, IH| , affect the efficiency of the 
EET. We also take into account the effect of the quan¬ 
tum nature of the environment, such as Markovianity 
and non-Markovianity on the dynamics [3. 135l - l37l| . 
Comparing the result of numerical calculation with re¬ 
cent experimental results, we analyze the effect of sur¬ 
rounding environment and inner spatial structure on the 
EET dynamics and efficiency. We find that the EET ef¬ 
ficiency is robust against various environmental parame- 
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ters. To shed more light on the quantum coherence that 
has been studied in many literatures in terms of entangle¬ 
ment [l3, [13, , we apply a recently developed mea¬ 

sure for quantum coherence called quantumness (43 - 1^ 
to study the quantumness of the EET within the single 
excitation subspace. By calculating the energy current 
of the transfer dynamics and its quantumness, we show 
that quantum coherence indeed play an important role 
in the energy transfer. 

The remainder of the paper is organized as follows. In 
Sec. II, we introduce a model to describe the complex 
quantum network, in which a phase factor is added into 
the coupling constant. This phase factor may come from 
the spatial structure of the light harvesting system. In 
Sec. Ill, we introduce a Markovian master equation in¬ 
cluding the relaxation, dephasing and dissipation to de¬ 
scribe the network. Then we apply the fully connected 
network to derive an exact expression for the EET effi¬ 
ciency. The dependence of the EET on the initial state 
and the phase factor in the model are discussed in de¬ 
tails. In Sec. IV, we focus on the Markovian and non- 
Mar kovian effects in the dynamics. The effects of the 
phase factor and temperature on the dynamics are also 
discussed. In Sec. V, we calculate the energy current in 
the system and discuss the quantum coherent nature of 
the dynamics. The effect of the phase factor and tem¬ 
perature on the current and quantumness are explored. 
Sec. VI is devoted to concluding remarks. 


II. MODEL 

We discuss a network of three sites, two of them |1) 
and 12) may support excitations, which can hop from 
one site to the other (see Fig. [T]). The two sites interact 
with the environment, and one of the sites (say site 2) 
is connected to the third site, sink |3). Once excitations 
fall into the sink, they can not escape. This network is a 
model simplified from the structure of the FMO complex, 
a network of seven coupled sites, each of which can be 
treated as a two-level system. The total Hamiltonian in 
the single exciton manifold includes three parts, H = 
Hs + Hsb + Hb- The system part is given by, 

Hs = eicr+ erf -I- £ 20 -^ erf -I- J erf -f e*‘^er^ erf) ,(1) 

where and erf are the raising and lowering opera¬ 
tors for site j, defined by cr+ = \j) (0| and erf = |0) (j|, 
where |0) denotes the zero exciton state of the system and 
\j) represents one excitation on site j. ej is the on-site 
energy of site j, and Jij is designated as the excitonic 
coupling between site i and j, where Jij = Jji = J is 
assumed. A phase factor with tp a real number is 
added to the inter-site coupling J, making the coupling 
constant complex. 

The Hamiltonian of the environment, which modulates 
site energy fluctuation and brings decoherence, is usually 



FIG. 1: (Color online) The schematic representation of three 
sites network. The exciton, initially in site |1), is transferred 
to site 12) and finally trapped by the sink |3). 

modeled as a set of quantum oscillators, 

Hb = ( 2 ) 

j m—1,2 

where b}^ j {bm,j) is the creation (annihilation) operator 
of the j-th mode with frequency ujj in the bath, we sup¬ 
pose that each site is separately coupled with its own 
environment. Thus the system-environment interection 
Hamiltonian reads [s^ 

Hsb = J2 97n,jibl,j + b„,,j)\m){m\, (3) 

j m—1,2 

where gmj denotes the coupling constant of the m-th site 
to its bath mode j. 

The transport dynamics of a single excitation is de¬ 
scribed by a master equation that includes coherent evo¬ 
lution, relaxation, dissipation and dephasing. In addi¬ 
tion, the exciton can be trapped in the sink. Once the 
excitation is trapped, it can not escape. 

III. EET DYNAMICS IN MARKOVIAN CASE 

Consider a Markovian environment, the dynamics of 
the reduced density matrix of the system is given by 

pit) = - i[Hs, p] + L^ipit)) 

+ Ljipit)) + Lripit)) + Lsipit)), 

where the Liouvilian on the right-hand side in Eq. 0 
describes coherent dynamics, relaxation, dephasing, dis¬ 
sipation, and sink, respectively. 

The relaxation (with rate ri 2 (T 2 i)) is described by 
Lindblad super-operators [13, [s^l , 

LK.ipit)) =ri 2 [cr^erf p(t)cr+erf 

1 + - VI + - 

-br2i[erfer2 /o(t)erferi 

- \ {pW,cr^cr)*-erferf }], 

where {A, H} denotes an anticommutator, while ri 2 de¬ 
notes exciting process from |1) to |2) and r 2 i denotes 
decay process from |2) to |1). 
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The dephasing process (with rate -jj) that destroys the 
coherence can be described by, 


same on every site. With those notations, the transfer 
efficiency is given by, 


J = 1.2 

Here, we treat pure dephasing in a way equivalent to the 
well-known Haken-Strobl model, in which pure dephasing 
is considered phenomenologically in terms of a classical 
fluctuating field MM- 

Additionally, two processes that lead to irreversible 
loss of energy are the dissipation and the sink. These 
effects are characterized by Lindblad operators [l3, [s^ 

Lr{p{t)) = Y ~ \ ^7}] (7) 

i=i,2 

and 

Ls{p{t)) =Ts[(T^<72 P(.t)<^t<^3 

Ir ... + _ + ( 8 ) 

- 2 {PW,<7^<73 aja2 |J, 

where Tj and Tg denote the rates of dissipation and trap¬ 
ping, respectively. In Ref. [H, [H, 113, , these two pro¬ 

cesses are taken into account by adding a non-Hermitian 
Hamiltonians Hioss = —ij^l X) |j) {j\ to the system 

i=1.2.s 

Hamiltonian o, which is equivalent to our treatment 
when the sandwich term can be ignored. 

The efficiency P of EET is quantified by the population 
transferred to the sink |3) from the site |2) 

P = P33{T)=T,[ P22{t)dt, (9) 

Jo 

which will be used as a measure of the transport effi¬ 
ciency. Thus, the question we focus on is as follows: In 
a given time T, what is the probability transferred from 
site |1) to the sink |3) and how the spatial structure of 
the system and the surrounding environment affect this 
transfer. 

To study how the coherent evolution and the noise- 
induced decoherence processes work together to reach ef¬ 
ficient EET and illustrate the role of the phase factor 
in the coupling constant, we adopt a uniform and non- 
uniform two-level ECN to help us identify some mech¬ 
anisms that underlie the energy transport [H, [l^. A 
uniform ECN means, ei =£2 =£, Ti = r 2 = r and 
7 i = 72 = 7- Given that the exciton is initially in site 
|1), then an exact analytical solution of P (see the Ap¬ 
pendix) for the case of T = 00 can be obtained 

Tg 

^ , -r \ , ^r(r+r,+r2i-ri2) ’ (10) 

-I- i gj -I- (4J2+^ri2) 


p= - 

( 2 r + rg) 


Tg 

r(A2-i-i6A2)(r-i-r„-i-r2i-ri2) ’ 
(4J^A-|-ri2(A^-|-16A^)) 


( 11 ) 


where the on-site energy gap is defined by A = 
(£2 — £i)/ 2. It is notable that, although the phase factor 
(f exists in the system Hamiltonian of Eq. CD, the final 
efficiency P does not depend on tp, which is in agreement 
with the findings in Ref. [l^. This is not the case, how¬ 
ever, if the exciton is initially in a superposition of | 1 ) 
and 12 ), as we show below. 

The dephasing process (represented by the term with 
7 j) does not assist the EET as expected, see Eq. (fTOj) . It 
is not the case, however, when we consider a non-uniform 
network, see Eq. CID, in this case P is not a monotonic 
function of 7 . We also observe from Fig. [5] (a) that, 
dephasing facilitates the transition only when A is larger 
than a certain value. The physics behind this difference 
is as follows. The range of on-site energies of |1) and 
12 ) are broadened owing to the dephasing, which leads to 
the overlap of the two sites in energy. If the dephasing 
rate continues to increase, the effect of resonant mode 
decreases as the energy of each site is distributed over a 
very large interval. When A is quite small, or even A = 0 
(uniform case), resonant mode is available without the 
presence of broadened site energy. 

From Eq. CHD or Eq. CID, it is easy to conclude that 
the dissipation rate Tj simply decreases the efficiency P 
(see Fig. [2] (b)). The dissipation process leads to irre¬ 
versible loss of energy via site j to the environment which 
obviously against the increasing of the EET efficiency. 

The influence of dephasing and dissipation in 
absence of relaxation (i.e.,ri2 = r2i = 0 ) on the 
transfer efficiency in a similar network has been 
exp lored extensively in the previous works [ii- 
IitL [ 2 ^. In the following, we are interested in 
new results obtained when the relaxation pro ces s 
(us ual ly ignored in Haken-Strobl model) (0- 
Il7l l22l| i is taken into account. We And from Fig.[^ 
(c) that P monotonically increases with ri 2 , but 
decreases with r 2 i, this finding can also be ob¬ 
tained from Eq. (IlOp or Eq. (Ilip by taking the 
derivative of P with respect to ri 2 or r 2 i. The 
observations can be understood as follows. Re- 
laLxation process represents thermal equilibration 
of the exciton. The diffusion of the exciton from 
|1) to |2) obviously increases the populat ion of site 
|2). If J = 0 is set, the efficiency P = 0 |l4l| . since 
no process except coherent transfer (with cou¬ 
pling constant J) transfers excitons from |1) to 
|2). When consider the relaxation process, the 
situation changes. For J = 0, Eq. (|10p becomes 


P = 


_TiaT^_ 

-I- r2irs -i- r (ri2 -i- r 2 i -i- r^) 


( 12 ) 


where A = 27 -|- 2r -I- Tg -|- ri 2 -I- r 2 i. For the case of This suggests that excitons may be transferred to 
a non-uniform ECN, £1 ^ £2 while 7 ^ and Tj are the the sink, even though coherent transfer character- 
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ized by J is zero. As the excitons can be trans¬ 
ferred via relaxation process, coherent transport 
and rel2Lxation process work together to get a 
highly EET efficiency. 

Fig- m (d) shows the dependence of the trans¬ 
fer efficiency P on the trapping rate and how 
these two processes work together to reach effi¬ 
cient EET. In the dashed, dotted and solid lines, 
the relaxation rates are the same, so the improve¬ 
ment in the efficiency is attributed to the coher¬ 
ent transfer. The red-dot-dashed and black-solid 
lines share the same values of J, but the relax¬ 
ation rates are different. The improvement of 
the efficiency is attributed to the relaxation pro¬ 
cess. Note that the thick orange-dashed and thin 
orange-dashed curves rise monotonically. which 
were not predicted in previous works |l4l- ll7l| . 
The other curves ascend first and descend then. 
This is due to the relaxation ignored in the lit¬ 
erature. This can be explained by analyzing Eq. 
(|10p . Taking derivative of P with respect to Fg 
(see Apendix B), we observe that, for a fixed F, 
ri 2 , F 21 and J, when < Fi 2 (F -I- F 12 -f F 21 ), P is 
a monotonic function of Fg, it takes a maiximal 
value P — with Fg —> 00, this is shown by 

the thick orange-dashed and thin orange-dashed 
curves in Fig. (d). In this case, the relarxation 
dominates the transfer. When the trapping rate 
Fg is very small, site |2) couples weakly to the sink 
|3), excitations can rarely reach the sink in this 
situation, so the efficiency is very low. When Fg 
is very large, the efficiency can be improved sig¬ 
nificantly. The highest efficiency is determined 
by the contribution of incoherence process given 
by P = When J2 > Fi 2 (F + F 12 + F 21 ), 

P is not a monotonic function of Fg and it also 
takes P = with Fg — 00, but there exists 

an optimized value of Fg (Fig. (d)), leading to 
a efficiency larger than From Fig. (d), 

the improvement of efficiency comes from coher¬ 
ent transfer. In this case, the coherent trans¬ 
fer dominates the transfer. When Fg is very 
small, the efficiency is also very low. When Fg 
increases, P increases gradually to a maximum 
where the contribution from coherent transfer is 
brought into full play and excitation transferred 
from |1) to |2) matches the transition from |2) to 
13) perfectly. Further increasing of Fg deterio¬ 
rates the coherent transfer, and eventually the 
efficiency decreases to the rel2Lxation-dominated 
case, i.e., P = when F, — 00. Thus, exci- 

tation transferred from |1) to |2) mismatches the 
transition from |2) to |3). In Ref |l4l| . only the 
latter case, i.e., > Fi 2 (F-|-Fi 2 -l-F 2 i), is discussed 

because the relaxation process is not considered, 
i.e., Fi 2 = F 21 = 0. 

The enhancement achieved by the dephasing can be 
understood as the fluctuation-induced-broadening of en¬ 


ergy levels, and this gives rise to a conjecture that the 
fluctuations in the site energies and couplings may affect 
the transfer efficiency. We will focus on this question in 
the following discussions. 

We add two energy fluctuations Ji, ^2, to the on-site 
energy gap A and inter-site coupling J separately. The 
expression of P Eq. m can be easily derived, 

_Ef_ 

/'ov I v ^ I r(A^-t-i6(A-t-(5i)^){r-t-rs) ’ (t-l) 

-I- i gj i- (4(j+52)^^+r,,(A2-|-16(A-|-5i)^)) 

where A = 27 -|- 2F -|- 2Fs -|- 2 Fk and F12 = F21 = F„ are 
performed. Obviously, it monotonically increases with 
1^2(or J), whereas it decreases with i5i(or A). This is easy 
to be understood. The energy gap A blocks the energy 
transfer, while the inter-site coupling J that represents 
the overlap of sites |1) and |2) favors the transport. Fig. [3] 
(a) shows this result. 

Now we examine the effect of phase (p on the 
transfer efficiency, ip in the coupling is determined 
by the spatial distribution of |1) and |2), this suggests to 
consider a superposition of |1) and |2) as the initial state. 
For the FMO complex, the spatial and temporal 
relaxation of excito n s hows that site 1 and 6 were 
populated initially [29| in sunlight. In laboratory. 
Ref reported the initial state preparation of 
the system in a coherent superposition of the an¬ 
tenna proteins electronic vibrational eigenstates 
with the femtosecond laser pulse (25-fs duration). 
Then it is reasonable to study how the initial states af¬ 
fect the excitation transfer efficiency. In the following, 
we shall shed light on this question and analyze whether 
the initial states have effects on the EET efficiency and 
whether the phase ip can enhance the transfer. We con¬ 
sider two different initial states. One is a superposition 
state of sites |1) and |2), and the other is a classical mix¬ 
ture of |1) and |2). 

We illustrate numerical results of the transfer efficiency 
as a function of a and /3, which characterize the pure 
initial state of the system through |-!/)(t = 0)) = cos a |1)-|- 
sinaexp(i/3) |2), where a denotes the population ratio 
between sites |1) and |2), and /3 characterizes the relative 
phase. From Fig. |3] (c), we observe that for excitation 
initially excited in a superposition of sites |1) and |2), 
both the population ratio and the relative phase affect the 
energy transfer, and suitable values of a and /3 facilitate 
the efficiency. The efficiency P reaches its maximal value 
when a approximately equals to ^ and This tells 
that occupation on site |2) helps the EET efficiency, see 
Eq. ([3]). Fig. |3] (d) shows the dependence of P on tp 
and a, and Fig. [3] (e) shows the dependence of P on 
and p. These two figures indicate that the phase ip does 
have effect on the efficiency when the initial state is a 
superposition of the two sites. This observation can 
be understood by an analytical solution 

P = 

Fg(8J^ -I- 2F12A -I- FAsi+a — 2JF sin2asin(/3 — ip)),^"^) 
(4J2(2F + Fg) + AF(F + Fg + F21 - F12)) 
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FIG. 2: (Color online) The transfer efficiency P. (a) P vs 7 and A. (b) P vs F. (c) P vs ri 2 and r 2 i. (d) P vs Fs. Unless 
otherwise noted, the other parameters chosen are F = 1, 7 = 0.1, J = 1, Fs = 1, ri 2 = r 2 i = 1, A = 0. In (d), the parameters 
is as follows: J = 0, ri 2 = r 2 i = 1 for thick orange-dashed, J = 0.5, ri 2 = r 2 i = 1 for thin orange-dashed, J = 1.5, 
ri 2 = r 2 i = 1 for thick green-dotted, J = 3, ri 2 = r 2 i = 1 for thin green-dotted, J = 4, ri 2 = r 2 i = 1 for black-solid, J = 4, 
ri 2 = r 2 i = 1.2 for thick red-dot-dashed, ri 2 = r 2 i = 1.4 for thin red-dot-dashed. 



FIG. 3: (Color online) The transfer efficiency P. (a) P vs and 52- (b) P vs p. (c) P vs a and /3 {p = ^). (d) P vs a and 
</9 (/? = j). (e) P vs P and p {a = !■). Unless otherwise addressed, the other parameters chosen are F = 1, 7 = 0.1, J — 1, 
Fs = 1, ri 2 = r 2 i = 1, A = 0. 
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obtained from Eq. Q with A = 0 and A = 27 + 
2r + Fs + ri 2 + r 2 i. The factor sin(/3 — ip) in the 
numerator indicates that varying the values of ip 
is equivalent to varying the values of the relative 
phase j3. This is actually equivalent to choose 
different coherent superposition of sites | 1 ) and 
12 ) as the initial states, for example, 

,, ( cos^a sinacosae”®^ 

p{t — 0 ) — I . .2 

ysinacosae '^ sm a 



Therefore, when the network is initially in a su¬ 
perposition state, the effect of tp on the EET is 
obvious. When p takes a fixed value, the rela¬ 
tive phase /3 also affects the coherence of exciton 
transfer and the transfer efficiency P. But when 
the exciton is initially in site |1), neither p nor /3 
has effect. 

In Ref. |47j |. the concepts of topological phases 
is introduced to explore the exciton transport in 
porphyrin thin films. This is the first work that 
addresses the joint effects of both magnetic field 
and coherence in molecular exciton transport. 
The magnetic field serves to break time-reversal 
symmetry and results in lattice fluxes that mimic 
the Aharonov-Bohm phase acquired by electrons. 
The first blueprint for realizing topological phases 
of matter in molecular aggregates stimulates us to 
study the EET dynamics from new point of view. 

Another type of initial states we consider is mixed state 
with mixing ratep, p |1) (1| + (1 — p) |2) (2|, see Fig.|3](b). 
We observe that p works like a denoting the population 
weight. More population on site |2) favors the transport, 
and p has no effect on the efficiency. This is because for 
sites | 1 ) and | 2 ) being initially in a classical mixed state, 
there is no coherence initially, leading to no influence of p 
on the efficiency. Actually, we obtain an exact analytical 
expression for P in this case. 


P= - 

(2r + r«) 


r. 

Ar[(2p-i)r+pr.] ’ 
A[(i-p)r+r„]+4j2 


(16) 


where A = 27 + 2 r + 2 rs + 2Tk. and ri 2 = r 2 i = F^ are 
performed. From this expression we easily know that P 
is independent of p. 

The effect of phase p on the EET dynamics and ef¬ 
ficiency will be discussed further in section V, using a 
time-convolutionless master equation. 


IV. EET DYNAMICS IN NON-MARKOVIAN 
CASE 


time-convolutionless master equation within the weak- 
cou pling, B orn-Markovian and secular approximations 
ksLl^. In contrast to the Haken-Strobl model, the 
environment is treated in a quantum way, 

p{t) = -i[Hs, p] + Lr{p{t)) + Lsip{t)) + L{p(t)), (17) 

where the first three terms on the right side are the same 
as those in Eq. ®, representing coherent evolution, dis¬ 
sipation and sink, respectively. L{p(t)) is the Lindblad 
superoperator given by [H, : 

Hpit)) = 

r 1 

Amiuj)p{t)Al^{uj) - - {Al,{uj)A^{uj),p{t)} , 


where is the non-Markovian decoherence rate un¬ 

der the assumption of Ohmic spectral density with expo¬ 
nential cutoff 

J(w) = —wexp ^ , (19) 

ujc \ / 


where lUc is the cutoff frequency and A is the strength of 
the system-bath coupling. Generally, the non-Markovian 
decoherence rate is time-dependent given by 


pCO 

7 (a;,t) = 2 / dQJ (ui) 
Jo 


, , sin [(w-I-w) t] , , sin [(w — w) d 

n (w)-+ [n (w) -f 1 ]- 


OJ + OJ 


( 20 ) 


where n (w) = [exp {uj /ksT) — 1 ] ^ is the bosonic distri¬ 
bution. In the Markovian limit {t —>■ 00 ), the decoherence 
rate can be obtained 


7 (w, 00 ) = 27rJ (|a;|) |n (—a;)|. ( 21 ) 

The dephasing rate can be derived from Eq. (I^Oll in the 
limit w —> 0 

Iv ^ { 2 ^') (22) 


In the Markovian limit {t —^ 00 ), the dephasing rate be¬ 
comes 


, . 2TrkBTX 

Iv ( 00 ) =- 

OJc 


(23) 


In the later simulation, unless otherwise noticed, we set 
A = 50 cm“^, OJc = 50 cm“^ as in Ref. [s^, and note 
that in units of /i=l, we have 1 ps“^ = 5.3 cm“^. 

Am{tj) is the jump operators defined by 


In section III, we treat pure dephasing in a way equiv¬ 
alent to the Haken-Strobl model, in which the envi¬ 
ronment is modeled as a spatially uncorrelated classi¬ 
cal white noise. Now we consider the EET dynam¬ 
ics using the Lindblad master equation derived from a 


Am(w) — C„^(^)Cm (0 (24) 

The sum runs over all possible transitions in the sin¬ 
gle exciton manifold. In our model, the basis |Afe) = 
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X] Cm (fc) |to) {k = 1 , 2 ) is composed by the eigenbasis 

m—1,2 

of Hamiltonian o, Hs\Xk) = Afe |Afc). Simple algebra 
yields, 


Ai — + £2 + \l (£1 ~ £ 2 )^ + 4J^], 

A 2 = ■ 5 [£i + £2 ~ (£1 ~ £ 2 ) + 4J^]. 

And the corresponding eigenstates are, 

|Ai) = sin0 |1) + cos 9 exp{i(p) |2), 
IA 2 ) = COS0 |1) — sin0exp(i</5) |2). 

where the mixing angle is given by 


(25) 


(26) 


sin 6 * = 


cosd = 


£1 ~ £2 + 4J2 + (ffi — £ 2 ) 



Thus the jump operators for relaxation are given by 


Ai ,2 (W 12 ) = ±sin6»cos6» |Ai) (A 2 I, 

Ai ,2 (w 2 i) = ±sin0cos0 jA 2 ) (All, ^ ^ 

for dephasing, they are 

Ai (0) = sin^e |Ai) (All + cos^6» IA 2 ) (A 2 I, 

A 2 (0) = cos^0 |Ai) (All + sin^0 jA2) (A2i. 


In the Hilbert space spanned by the eigenstates of the 
Hamiltonian, L{p(t)) can be written as 


7^ {p (0) — ^ ^ ^ni^nnP (f) ^nn ^ P (7)}) 

n=l,2 

+ 712 (ctuP (t) (722 + CT22P (i) CTii) 

+ r2i(cr2ip(<) CT12 - ■^Wll,pit)}) 

+ Ti2{ai2p{t) (721 - ^ {ct22,P (<)})> 

where (Jy = |Ai) (A^j. The first two terms represent the 
dephasing and the last two terms describe relaxation, 
which is similar to Eq. jS]) and Eq. The relaxation 
and dephasing here is expressed in the basis composed 
of the Hamiltonian eigenstate, different from the phe¬ 
nomenological treatment in Eq. ([5]) and Eq. ([ 6 ]). The 
decoherence rates defined in Eq. (uni-nsi) take, 


7 „ =(sin'^6< -I- cos'^ 6 ») 7 ,^ (t ), 
7 i 2 = 2 sin^ 0 cos^ 07 y,(t), 
ri2 = 2 sin^ 0 cos ^07 (a;i2,t) , 
r2i = 2 sin^ 0 cos ^07 {uj2i,t) i 


for non-Markovian dynamics. They take, 

7 „ =(sin'*0 -I- cos‘^9)j^ (oo), 

712 =2sin^0cos^6i7,^(oo), 
ri 2 = 2 sin^ 6 icos^ 6*7 (a;i 2 ,oo), 
r 2 i = 2 sin^ 0 cos ^07 (a; 2 i, oo), 

for Markovian case, where ojki = Xk — is the transition 
frequency. The parameters chosen are: J = 87 cm“^, 
£ 1=0 and £2 = 120 cm“^, equivalent to the Hamiltonian 
of sites |1) and |2) subsystem in the EMO complex in Ref. 

M- 

With these arrangements, we numerically calculated 
the time dependence of the decoherence rates, the time 
evolution of the population on the three sites and the 
coherence element of the density matrix in the site 
basis, i.e., Pmnit) = (to| p (t) |n). Markovian and 
non-Markovian decoherence rates and transfer dynamics 
given by Eqs. (UHl- HSl) are plotted in Eig. [Hand Eig. [S] 
The non-Markovian relaxation rates oscillate, taking pos¬ 
itive and negative values, and finally converge to values in 
the Markovian limit [s^ . The dephasing rates in the non- 
Markovian case start with zero and similarly, converge to 
that in the Markovian limit. We also observe that higher 
temperature leads to larger amplitude of oscillation in 
the relaxation rates. Eor the dynamics, the oscillations 
of population in sites | 1 ) and | 2 ) and coherence element 
in the Markovian case at low temperatures are almost the 
same with those in the non-Markovian case. At higher 
temperatures, however, considerable differences can be 
observed from Fig. [Hand Fig. [S] Oscillations in the pop¬ 
ulation and coherence element last shortly, whereas in 
the non-Markovian case, they last long. This can be 
explained as follows. For Markovian case, the dephas¬ 
ing rate is linearly proportional to the temperature ac¬ 
cording to Eq. (|23)) . whereas at the same temperature, 
non-Markovian dephasing rate is much smaller at short 
time from Fig. |4l In addition, the Markovian relaxation 
rate is always positive, while for the non-Markovian case, 
the relaxation rate may take positive and negative. Since 
the relaxation process represents the diffusion of excitons 
through the system, oscillatory relaxation rates of non- 
Markovian dynamics lead to beatings of the population 
and coherence element, which reflects the memory effect 
of the non-Markovian environment. 

These factors lead to the difference of the pop¬ 
ulation on site |3) between Markovian and non- 
Markovian cases. Note that we consider the pop¬ 
ulation P 33 accumulated on the site 3 for a time of 
t — 1 ps, i.e., we choose P 33 (1 ps) to quantify the 
transfer efficiency, this is different from the case 
in Sec. HI, where we use p 33 (oo). Similarly, at low 
temperatures, the population in site |3) in Marko¬ 
vian case is almost the same as in non-Markovian 
case. At higher temperatures, it exhibits a slight 
difference (see Fig. |4] and Fig. [ 6 ] (a)), which can 
also be explained by the different behavior of 
decoherence rates between Markovian and non- 











Markovian 


Non-Markovian 






Time[ps] 



FIG. 4: (Color online) Time-dependent decoherence rates from Eq. (I20I21I22I23II and time evolution of population on each 
site from Eq. G3 for Markovian and non-Markovian dynamics. The left panels (Fig. 0] (a) - (d)) show the time evolution of 
relaxation (upper panel) and dephasing (lower panel) rates. The parameters in the spectral density are Uc = 50 cm“^ and A = 
50 cm“'^ as used in Ref. |35||. The relaxation rates (blue-dotted for uJi 2 and green-dot-dashed for uj 2 i) oscillate, they may take 
negative value and converge to values in the Markovian limit (red-solid for 0 J 12 , the purple-dashed line is for UJ 21 ). Similarly, 
the dephasing rates in the non-Markovian case (blue solid) start with zero and converge to these in the Markovian limit (red- 
dashed). Fig. 0] (e) - (1) show the population on each site as a function of time: poo in purple-dashed, pn in blue-dotted, 
P 22 in green-dot-dashed, and P 33 in red-solid, (e) - (h) is plotted for the Markovian case, and (i) - (1) for the non-Markovian 
case. The temperatures T and phase ip change from figure to hgure. The population on the reaction center P 33 with the 
same temperatures T is also illustrated (red-solid for tp = 0, orange-dashed for ip = 7 r/ 2 , black dot-dashed for p = tt). The 
Hamiltonian is taken to be the same as that in Ref. [l^ . Dissipation and trapping rates we chose are Fi = 0.1 ps“^, r 2 = 0.1 
ps“^, and Fs = 10 ps“^. 
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FIG. 5: (Color online) Time evolution of the coherence elements (Re pi 2 red-dashed, Im pi 2 blue-solid) from Eg. (1171) for 
Markovian and non-Markovian dynamics. The excitonic Hamiltonian is taken to be the form equal to that in Ref. |13l] . The 
parameters are A = 50 cm~^, u)c = 50 cm“^, phase ip = 0, and room temperature T = 300 K. Dissipation and trapping rates 
we choose are Fi = 0.1 ps“^, r 2 = 0.1 ps“^, and Ts = 10 ps“^. 



Temperature [K] 

(a) 



Reorganization energy A. [1/cm] 

(c) 



(b) 



(d) 


FIG. 6: (Golor online) The transfer efficiency P as a function of different temperatures, the phases p and the main decoherence 
parameters (reorganization energy A, cutoff uJc) for Markovian (red-dashed) and non-Markovian (blue-soloid) cases. The 
parameters are ujc = 30 cm“^, A = 30 cm”*^, phase = 0, which is typical for some natural energy transfer systems like 
chromophores in photosynthetic systems, and the temperature is T = 300 K. Dissipation rate ri, 2 , trapping rate Fs, and 
excitonic Hamiltonian takes the same values as in Fig. [d] In (b), thick red-solid for Markovian and thin blue-solid for non- 
Markovian with the initial state pii(O) = 1 , and red-dashed for Markovian and blue-dot-dashed for non-Markovian with the 
initial state a superposition state \tp{t = 0)) = | |1) + -^ |2). 
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Markovian cases. Note that the difTerence in the 
efficiency between Markovian and non-Markovian 
cases is not so large. Furthermore, from Fig. [6] 

(a) , we find that the temperature does not change 
much the population on site |3), hence the effi¬ 
ciency is robust against the temperature. This is 
consistent with recent experimental results that 
the nature photosynthesis are able to transport 
incident light energy with n earl y 100% efficiency 
in hot and wet conditions f3-[Tol|. 

Fig. [6] (b) shows the efficiency as a function of 
the phase ip for Markovian and non-Markovian 
cases. Fig. [6] (b) show that the phase factor ip 
has influence on neither the EET dynamics nor 
the efficiency with initial state /Oii(O) = 1, this 
observation holds true for both Markovian and 
non-Markovian cases. When we consider a su¬ 
perposition of |1) and |2) as the initial state, the 
effect of ip on the efficiency is clear from Fig. [6] 

(b) , and the improvement in the efficiency com¬ 
pared with pii(O) = 1 is obvious. These results are 
in agreement with the conclusions in section III. 
The effect of the phase factor on the EET depends 
on initial states. The difference in the efficiency 
for Markovian and non-Markovian cases is very 
small. 

In Fig. [5] (c), we show the dependence of transfer ef¬ 
ficiency on the reorganization energy A of the spectral 
density. From Eq. (uni-iMi), we see that relaxation and 
dephasing rates are both proportional to A. As discussed 
above, the dephasing destroys the coherence and hinders 
the energy transfer when the on-site energy gap A, i.e., 
Wi 2 (w 2 i) in this section, is sufficiently small compared 
with the dephasing rates, while it improves the energy 
transfer when A is suitably large compared with the de¬ 
phasing rates. But the relaxation terms describing the 
exciting and decay processes among the eigenstates |Ai) 
and IA 2 ) would only accelerate the energy transfer. With 
A getting larger from a very small value, both the dephas¬ 
ing rates and relaxation rates become larger. Therefore, 
the efficiency varies nonmonotonically with A increasing 
for both Markovian and non-Markovian dynamics. The 
difference in the EET dynamics for Markovian 
and non-Markovian cases results from the differ¬ 
ent decoherence rates when A changes (Fig. [4] (a) - 
(d)). The difference in the transfer efficiency can 
be even small, see Fig. [6] (c), this is because the 
transfer efficiency is defined as an accumulation of 
population on site 3. It is worth addressing that 
the master equation (ED is obtained within the 
weak-coupling limit. Therefore, for small values 
of A, the numerical results is reasonable. When A 
is large, the master equation (|17|) is theoretically 
not applicable, this means it can not be derived 
as in textbook, but it is allowed to use the master 
equation phenomenonaly to describe the EET. 

Fig. [S] (d) shows the dependence of P as a function 
of (jJc- Eq. (0 tells that the spectral density J{w) is 


not a monotonic function of Wc, and the four decoher¬ 
ence rates change differently with lOc- Similarly, when 
LOc changes, the roles that the relaxation and dephas¬ 
ing play in the dynamics change for both Markovian 
and non-Markovian dynamics, leading to the difference 
in the dependence of P on ojc between Markovian and 
non-Markovian cases, see in Fig. |n](d). Notice that a 
larger cutoff ujc means stronger couplings between 
the system and environment, so the master equa¬ 
tion (jl7D becomes not applicable from the respect 
that it can be derived by the use of week-coupling 
approximation. 

The reason why the different surrounding environ¬ 
ments can affect the EET dynamics but exhibit no in¬ 
fluence on the efficiency (Fig. 2] - Fig. IHD lies in that: 
when the EET dynamics is in the stage of coherent os¬ 
cillation, the exciton is transferred relatively faster. The 
surrounding environment and the temperature may influ¬ 
ence the duration of coherent oscillation, but before this 
stage ends, the exciton has been transferred to the sink 
13) with a large probability, no matter how long this pe¬ 
riod of coherent oscillation is. This suggestes that quan¬ 
tum coherence play a significant role in achieving the re¬ 
markable efficient EET The explanation above 

gives an intuitive picture how the exciton is transferred 
through the system which will be further illustrated in 
terms of the notion of quantumness in Sec. V. 

It is notable that the efficiency of EET for Markovian 
and non-Markovian cases discussed in exhibits con¬ 
siderable difference. It is due to the different measure 
to elucidate the EET efficiency. In this paper, we adopt 
trapping sites |3) (with trapping rate Eg) to model the 
sink while [s^ ignore the trapping process and utilize 
another kind of measure given by 

P=- f dt{M\p{t)\M), (33) 

^ Jo 

where t is the total integration time and \M) is a partic¬ 
ular site. 

V. QUANTUMNESS OF ENERGY CURRENT 
IN EET 

In Sec. IV, we have shown that quantum coherent 
wave-like oscillations of the populations in sites |1) and 
|2) last up to 500 fs at temperature T = 10 K, and be¬ 
come shorter with the temperature increasing for both 
Markovian and non-Markovian dynamics. In the follow¬ 
ing, we will quantify how quantum the EET dynamics is 
and study the influence of environment and spatial struc¬ 
ture of the system on this process. As energy transfer 
constitutes the main function of the system, we focus on 
the quantumness of the energy transfer current . As is 
widely accepted, the most classical states of a quantum 
system are the pointer states, which are einselected by 
the decoherence process. For the case of a measurement 
of energy current, the pointer states are the eigenstates 
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of the energy current operator [4^, |4^. We then can 
measure the distance of the density matrix to the convex 
set of the pointer states [ 13 , [HI, i-e., classical mix¬ 
ing the pointer states, and define the quantumness as the 
minimum distance. 


A. Energy currents in the transport 


Comparing with Eq. dSzl), we obtain the energy current 


Sl->.2 = 2 ^ (^12^2/ |1) (/| — H.C.), 
i 

S 2^1 = ^^(h2ihii|2) {l\-H.c.). 

I 

The energy current operator is 


(40) 


First, we derive the energy current operator for our 
system utilizing a well-developed theory by Hardy 
which was wid^ used by a large amount of succeeding 
literatures [54l458l| . An energy current operator j (x) can 
be obtained from a continuity equation 

i7(x)-PV-r(x) = 0, (34) 

with H (x) the local energy density. For the Hamiltonian 

2 

m in tight-binding approximation, H = X) ^mn \iti) (n-|, 

m,n 

we rewrite the decomposition of the Hamiltonian in terms 
of local excitations as 

H — hi -\- /i2, 

1 A (35) 

= 2 ^ (^mn I™) {A + H.C.), 
n—1 


Consider a one-dimensional chain, the current can flow 
only in one direction, hence j has only one component, 
which will be denoted by j, and its value on site m de¬ 
noted by jm- Taking as equidistant with lattice constant 
a, we obtain the discretized form of Eq. (l34l) 

r) h T* — t’’ 

"'rn _ Jm Jm J’ 30 ^ 

dt a a ’ 

with /m^ the energy flux at positive x direction on the 
left (right) of site m. The balance of the currents si ^2 
from site |1) to site |2) and S 2 ->.i from site |2) to site |1) 
forms current jm ■ With these observations, we have 


— hi — 52^1 — Sl^.2, 
at 


— ^2 — Sl^2 — S2^.1- 

at 


(37) 


At the same time, we calculate the left-side of Eq. dSIl) 
by means of Heisenberg’s equation of motion, 

^^hm=l[H,hm]. (38) 

We have 

^ E |2) (^1 - hi2h2l |1) {l\ - H.C.), 

n ■ ' (39) 

(*12/^2/ |1) (^1 - h^ihii |2) {l\ - H.C.). 


jl 2 = Sl _>2 — 52^1 

~ 2(^11 ^22)(^12 |1) (2| — ^21 |2) (1|). 

Therefore, 

(ji2(t)) =tr[ji2p{t)] 

= 2 (ei + £2) J (e *‘^p2i {t) — e*'^pi2 (t)) , 


(41) 


(42) 


where pi 2 {t) and p 2 i it) are the time evolution of the 
coherence elements obtained from Eq. (El). We numeri¬ 
cally calculate (ji 2 [t)) at different temperatures for both 
Markovian and non-Markovian dynamics, and present 
the results in Fig. [71 The result of (ji 2 [t)) is given in 
unit of 6.57 x 10“^ eV/ps. As phase p does not affect 
the EET dynamics given the initial state pii(O) = 1, we 
don’t consider the influence of p. 

Lower temperature leads to longer oscillations and 
larger oscillation amplitude for both Markovian and non- 
Markovian cases, and the amplitude difference between 
the two cases become more obvious at larger tempera¬ 
ture. This is because the coherence die out very fast 
and more differences between the Markovian and non- 
Markovian dynamics emerge at high temperature, which, 
from Eq. (H7|) . will lead to smaller oscillations and larger 
amplitude difference. After the initial time, energy cur¬ 
rent becomes small but finite up to 1 ps, which helps in 
understanding the EET process further based on the con¬ 
clusion obtained in Sec. IV as follows: energy transfer 
begins with quantum coherent population exchange be¬ 
tween sites |1), 12) and |3), then small currents gradually 
redistribute energy between the three sites and finally the 
system arrives at an equilibrium. 


B. Quantumness 


The eigenvectors of ji 2 , i-e., the pointer states 
are given by 



defining the relevant pure classical states for the energy 
transfer. The minimum distance of p (t) to the convex 
set of classical states is defined as quantumness of the 
current ji 2 [soj . 


Qi 2 ip)= ^ , min 

{pi\pi>0,J2iPi=^ } 




Z=l,2 


,(44) 
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FIG. 7: (Color online) Energy current {ji 2 {t)) (top) and corresponding quantumness Qi 2 {p) (bottom) for the case of Markovian 
(red-solid for (ji 2 {t)) and blue-dashed for Qi 2 (p)) and non-Markovian (black-dashed for {ji 2 {t)) and green-solid for Qi 2 {p)) 
dynamics with different values of temperature T. Dissipation rate ri, 2 , trapping rate Fs, spectral density parameters A, lJc, in 
the excitonic Hamiltonian take the same values as in Fig. |4] 


where ||a 1|| is the Hilbert-Schmidt distance ||a1|| = 

(trAA^^ . Note that this measure of quantum¬ 
ness is introduced to characterize the bipartite 
quantum correlations between spatially separated 
sites |1) and |2), not the c orre lation between the 
system and environment |24l |44| . l50l| . The ef¬ 
fect of the environment on the quantumness is 
included in the density matrix p given by Eq. 
(|17p . According to this definition, Q 12 (p) > 0 where 
Qi 2 (p) = 0 if p is classical. An upper bound is given by 
Qi 2 (p) < Qmax = \/trp‘^ — 1/d, with d the dimension of 
Hilbert space [3. In our system, Qmax — 0.707. A finite 
value of Qi 2 (p) indicates that there are coherences left 
in p (t) written in the pointer basis. We numerically cal¬ 
culated the quantumness Q 12 (p) for the current {ji 2 (t)) 
at different temperatures for both Markovian and non- 
Markovian dynamics, and present the results in Fig. [3 
where Q 12 (p) is in unit of 1. Either, we don’t consider 
the influence of (p, as we choose pii(O) = 1 as the initial 
state. 

Similarly, the quantumness Q 12 (p) shows oscillations 
at the beginning of evolution, especially for lower tem¬ 
perature and non-Markovian case, but it decays fast at 


higher temperature. This also can be understood as 
rapid disappearances of the coherent oscillations. It is 
notable that Q 12 (p) decreases much slower than the os¬ 
cillation of the energy current. From Fig. [3 for T = 10 
K, Qi 2 (p) drops from its initial value (Qmax — 0.707) 
to Qi 2 — 0.3 — 0.4 within ^ 200 fs, but slowly drops to 
Qi 2 ~ 0.1 at ~ 800 fs. From Fig. S] and Fig. [3 the pop¬ 
ulation and the off-diagonal elements of the density ma¬ 
trix does not exhibit coherent features at ~ 800 fs, while 
the superposition of the pointer states of the energy cur¬ 
rent operator retains substantially coherent, comparing 
with the quantumness of the thermal equilibrium state 
Pt = exp{—Hs/kBT) at the four values of temperature 
adopted above, which all go to zero as the irreversible 
loss of energy to the sink and environment. That is to 
say, energy or exciton transfers through the network to 
a large extend in the form of a coherent snperposition of 
pointer states of the energy current operator. Quantum 
coherence plays a significant role in achieving the highly 
efficient EFT. [uHil- 
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VI. CONCLUSION 

In conclusion, we have studied the quantum dynamics 
of EET through a complex quantum network. We show 
how the coherent evolution and the environment-induced 
decoherence, especially the relaxation, work together for 
the efficient EET. The effect of the phase factor in the 
coupling constant on the EET depends on initial states, 
e.g., when the exciton is initially in site | 1 ), the phase 
factor ip has no effect on the transfer efficiency, but when 
the initial state is a superposition of | 1 ) and \2), ip really 
has effect on the efficiency. Using the Lindblad mas¬ 
ter equation, we evaluate the time evolution of popula¬ 
tion on each site for both Markovian and non-Markovian 
cases, we find that the dynamics depends on the Marko- 
vianity of the system but the efficiency is robust against 
the parameter variations of the environment. Finally, 
we quantified the quantum nature of the EET dynamics 
by calculating the energy current in the network and its 
quantumness, which helps understand the EET process 
further. It has been found that the energy current man¬ 
ifests substantial quantumness in both Markovian and 
non-Markovian dynamics. So, to some extend, we can 
say the EET dynamics is coherent despite its coupling to 
environments. 
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Appendix A: The derivation of Eq. (113 

In terms of the density matrix elements in the site basis 
Pij (t ), the equation of motion are 

pii(t)=-f(Je-*V 2 i-Je*Vi 2 ) 

— Tipii — ri2Pii + T21P22, 

P22{t) = — — Je“*‘^p 2 i) 

— r2P22 + ri2Pll — r2lP22 — rs/322, 

P33(i) =rsP 22 , 

poo(i) =ripii -I- r 2 P 22 , 

Pl2{t) = — *[—2Api2 -I- Je“*‘^(p 22 — Pll)] 

— 2 + r 2 -f 7i -I- 72 + ri 2 -I- r 2 i -i- rs)pi 2 , 

where A = (^2 — £i)/ 2 , and the initial conditions are 

Pii(O) = l,poo(0) = P22(0) = P33(0) = Pl2(0) = 0.(A2) 


By means of Laplace transform, the coupled differen¬ 
tial equations can be converted into a set of algebraic 
equations for the Laplace s-domain variables, 

shi=-i{Je-^^P2i-Je^'^Pi2) 

— Lipii — ri2pii -I- r2ip22 +1, 
sp22 = -i[Je^^Pi2-Je-^^P2i) 

— T2P22 + ri2Pll — r2lP22 — ^sP22, 

SP33 =rsP22, 

spoo =riPll -I- r2P22, 

spi2 = — *[—2Api2 -I- Je“*‘^(p22 — Pll)] 

— 2^^! + r2 -I- 7 i -I-72 + ri2 -I-r2i -i-rs)pi2. 

From Eq. (R3l) . we can easily obtain the expression of 
P 33 {s), and the relation of the Laplace transform for s 
and t gives 

p = P33ioo) = lim sp33{s). (A4) 

s->0 ' ' 

Then we obtain the expression of Eq. (fTTl) . Eq. (fTUl) . 
Eq. (fT^ . and Eq. (fTHll . Eq. (fT4)) can be got in the same 
way. 


Appendix B: Derivative of Eq. (I10|) with respect to 

r. 

Now we take the derivative of Eq. m, and obtain 

dP B + CTs+DTs^ 

df; =- E -' 

where 

C=2rri2(8j2 + (r + ri2 + r2i), 

D =r(-4j2 + ri 2 (r + ri 2 + r 2 i)), 

E = [{AJ^ + ATi2){2T + r,) (B2) 

+ Ar(r + r,-ri 2 + r 2 i)]", 

A =27 -i-2r -i-rs-i-ri2-i- r2i, 

and i? > 0 is a very complex expression without T^. 
From these expressions, the denominator E > 0. In the 
numerator, C > 0 and from the quadratic coefficient D = 
r(-4j2-pri2(r-tri2-br2i)), when .P < ri2(r-fri2-f 
r2i), D > 0, thus P is a monotonic function of Tg. when 
,P > ri2(r-|-ri2-l-r2i), U < O, then P is not amonotonic 
function of Tg, and it increases first and then decreases. 
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